This invention relates to methods for determining the relatively long time dynamics of molecular systems. The method also predicts the conformations of macromolecular systems including those that have undergone mutation.
Predicting how molecules move and what conformations they adopt is a problem that has important consequences in a variety of commercially important technical areas. For example, new drug development increasingly relies on the rapid prediction of molecular conformations to identify a few promising candidate compounds. By identifying from a large pool of candidate compounds those few possessing conformations consistent with a desired activity, the researcher saves considerable time that would be lost synthesizing and testing many different compounds.
One widely-used procedure for simulating the motions and conformations of molecules, and especially proteins, is molecular dynamics. Molecular dynamics computations are now extensively used for refining molecular structures obtained by X-ray or NMR techniques and for calculating the free-energy differences essential to correct evaluation of binding equilibria and the changes introduced by site-specific mutagenesis. Generally, in a molecular dynamics simulation, one expresses all the forces tending to change the positions of atoms within the molecule and then integrates Newton's force equation (F=ma) to obtain velocities and positions of the atoms in the molecule at a specified temperature. The relevant forces employed in a molecular dynamics simulation include both interatomic covalent and non-covalent forces. By repeatedly moving the atoms to new locations and then integrating over an appropriately short length of time (the "time step"), researchers can predict the motions of atoms within a molecule over a relatively short time frame.
The dynamics of a macromolecule comprises an enormous range of time scales including atomic vibrations on a subpicosecond time scale, amino acid rotational isomerization on a nanosecond time scale, and nucleation for helix formation on a 100 ns time scale. Molecular dynamics often provides a quite accurate description of the faster events (e.g. vibrations on a covalent bond) but can not be extended to the very much longer time scales where many events of commercial importance occur. This is because the accuracy of molecular dynamics simulation relies on careful numerical treatment at the level of the fastest time scale in the problem, thus requiring too much computational effort to reach the longer time regimes. This especially true when solvent molecules are explicitly included in the model, as is sometimes necessary to accurate describe motion. The time steps employed in molecular dynamics simulations are typically of the order of a few femtoseconds. Thus, given current computational power, a few 10 s of picoseconds (up to possibly 1 nanosecond) is the longest realistic time frame for molecular simulation. Unfortunately, this falls far short of the time domain over which such interesting large scale events as protein folding, opening fluctuations, and helix-coil transitions take place. Even with advances in computer technology such as progress in parallel computing, molecular dynamics may still be unable to describe such large-scale motions of macromolecules.
In order to overcome these limitations, it becomes necessary to introduce simplifying assumptions about the forces and the motions they cause within molecular systems. If these assumptions are formulated approximately correct, a sufficiently accurate simulation can be preserved while the simulation is extended into much longer time domains. One widely used assumption is that some of the fast modes of the molecule can be locally averaged and described as random noise forces that are balanced by frictional forces. This approach is gaining popularity and is especially important in simulations accounting for the effect of solvent on the molecular motions (i.e., Brownian motion). Because the noise forces introduce a statistical component, these simulations have been referred to as "stochastic dynamics." They are also often referred to as "Langevin dynamics" because the relevant force expressions containing frictional and noise forces are known as Langevin equations of motion.
Unfortunately, most stochastic dynamics methods still explicitly consider all the intramolecular motions, fast and slow, considerably limiting applicability of the simulation to long time molecular events. For example, Bhattacharya et. at. in International Journal of Quantum Chemistry, 42:1397-1408 (1992) describe a stochastic dynamics model of a small protein (bovine trypsin inhibitor) in which covalent and non-covalent forces are explicitly considered. The simulation employs a Langevin dynamics expression including covalent constraining forces and relatively complex frictional terms to account for the effects of water molecules. The authors claim that their method is one order of magnitude faster than a comparable molecular dynamics method, but this is still too slow to observe many interesting events.
Because the covalent forces associated with molecular systems are responsible for the faster motions such as bond vibrations, it might be useful to remove them from explicit consideration in the simulation. A 1977 paper (van Gunsteren and Berendsen, Mol. Phys., 34:1311-1327 (1977)) describes macromolecular dynamics method in which the covalent bond forces are replaced with constraint forces which simply maintain the relative positions of atoms bonded to one another within the molecular representation. In the molecular modeling community, a widely-used version of this procedure is known as "SHAKE." The method employs Newton's motion equations to model the motion resulting from the relevant forces. Like the other approaches to molecular and stochastic dynamics, these equations are second order (i.e. they include terms having second derivatives with respect to time). Thus, the method must determine not only the force acting on each atom but also the change in velocity. Such equations must then be solved by iteration (and therefore many matrix inversions) at each time step. This considerably complicates and slows the process.
Recently, Deutch et al. (Journal of Chemical Physics, 90:2476-2485 (1989)) described a stochastic technique for modeling the movement of DNA during gel electrophoresis. The authors modeled the DNA as a thread with beads and the gel as an infinite lattice of obstacles through which the thread and beads moved. They used a Langevin force description including terms for friction between the bead and a solvent, a constraining force between adjacent beads, a force between the beads and the obstacles, a random force acting on the beads, and a force of an applied electric field. This work is interesting in that it employs constraining forces in the Langevin description of molecular motion and employs first order equations which can be solved rapidly. However, it treats the DNA molecule as an extended chain of beads which do not interact with one another, other than through the constraining forces. In other words, it fails to account for the non-bonding interatomic potentials between the atoms comprising the DNA molecule. While this may be an adequate assumption for some macromolecular movements, for many other systems such as proteins and peptides, molecular movements and conformations are largely determined by non-bonding interatomic forces. Thus, the approach of Deutch et al. is inadequate for describing the dynamics of many important macromolecular systems.
Although conventional molecular dynamics and stochastic dynamics methods are useful in describing molecular events occurring on fast time scales, it would be desirable to have accurate methods of describing events on slower time scales.